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Abstract : 

This contribution describes the evolution of the protoplanetary disk using 2D numerical simulations. 
The 2D Euler equations are solved with the finite volume method. The numerical simulations are used 
to study the persistence and migration of anticyclonic vortices. Two cases are presented: (1) vortices 
produced by the Rossby Waves Instability, (2) a non-linear vortex model initially implemented into the 
disk. The migration of the vortices is due to spiral density waves excited by the vortex in the gas of 
the disk. 
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1 Introduction 

To understand protoplanetary disks' evolution is crucially important throughout the problematic of 
planetary formation. Therefore we are confronted with two major problems. The first appears to 
be a problem of time scale. Indeed the observations of protoplanetary disks and especially in the 
infrared excess resulting from dusts lead us to believe that the disk is dispersed after 10 7 years Q3J). 
The second problem is the growth of particles inside the disk. Actually several theories have been 
proposed to explain how some micron sized particles present a the beginning of the disk's life can 
reach a kilometer size and form planetesimals or planetary cores ([12], E] Gj The principal boundary 
to the grain's growth is that the flow is sub-Keplerian whereas the motion of particles is keplerian. 
When the size of solid particles is large in front of the mean free path of gas molecules, of the order of 
10 cm to 1 m, the particles undergo a wind which makes them fall onto the star in 100 years typically. 
The coagulation or particle sedimentation's time is longer by at least one order of magnitude. 

Hence the particle's confinement in an environment where this radial migration is reduced, is favorable 
for the particle's growth within a time scale sufficiently short to create gaseous giants. Lots of studies 
bear on anticyclonic vortices because they capture particles in their core. ([I], [5]) Cyclones being 
quickly dispersed by the Keplerian shear, only anticyclonic vortices do survive longer enough in the 
flow. Several processes make it possible to form anticyclones, and particularly the Rossby Waves 
Instability or the baroclinic instability which lead to the formation of large anticyclonic structures. 

We propose to study those vortices in a 2D model, based on the Minimum Mass Solar Nebula Model 
which is the standard model of protoplanetary disks. In a first part, we will clarify the cylindrical 
Euler equations, and then a quasi linear approach of the Euler equations allows us to propose a vortex 
model. In a second part, we will focus on different numerical simulations, with the study of several 
Rossby Wave Instability simulations and vortex models. To finish, we will explain those different 
simulation results, especially the vortices life time and their migration. 

2 The Euler equations and a vortex model 

The standard model of 2D protoplanetary disks is based on the study of matter repartition in the 
Solar System and on disks observations. We suppose that density, temperature and pressure decrease 
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like a power law of the radius. We assume that there is an hydrogen diatomic gas, polytropical with 
the following conditions in the thermodynamical equilibrium: 

- 1700 (ilt/)" s an ~ 2 T « - 280 (uu)' K 

P k B T a p + g = 7|) 
with — 2 < p < and — 1 < q < 0. We will choose here p = —3/2, and 7 = 1.4. 

2.1 The Euler equations 

The disk gas is in rotation in the Keplerian potential of the star of mass M*. The Euler conservation 
equations make it possible to describe this gas evolution. In the cylindrical coordinates, we write 
down, with V = u u r + v uq\ 

d t a + -d r (rau) + -d e (av) = (1) 

TV 
1 1 V 2 V 2 

d t au + -dJrau 2 ) + -d e (auv) = a g-± - d r P (2) 

r r r r 

d t crv + -d r (ravu) + -do(av 2 ) = —a— d$P (3) 



1 1 v 2 

d t cre + -d r (r(ae + P)u) + -d e ((ae + P)v) = -au^ (4) 



with Vk = y and the total gas energy ae = + |(^ 2 + v 2 ). 

A stationary solution associated with thermodynamical quantities ctq and Po is given by the following 
velocity Vo = uo u r + vo uq\ 

^o = (5) 
Po\ 1/2 



v = \ v k + — J (6) 

where the second term of the right member derives from ^rd r Po- We can notice that the gas is in 
sub-Keplerian rotation because of its pressure. We will call it the stationary solution. 

2.2 A vortex model 

We search a quasi stationary non-axisymetrical solution of the Euler equations for this flow. We set 
down: 

V = Vo + V (7) 

V = u u r + v uq (8) 

a = ao a (9) 

P = PoP (10) 

= £ = 4of (ii) 

If we get into a rotating frame with Vq, the continuity and velocity equations become: 



v 

dtcr + u d r a H — d$a + a 
r 



-d r (ru) H — d$v 
r r 



2^ 2 

+ ^ a r ^x H — o$u — v — h v c Sq 



r r r 



= (12) 

(13) 



r \ cr / a 



+ 14 a r t; + -a^ —11 ^ ^-deP (14) 



r r r a r 
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We search a stationary solution, so we lay dt = 0. If we suppose the terms in uv, v 2 , u d r and ^<9# to 
be small in front of the other terms, then the velocity field is given by: 



M: 



r 1 



u -v 



d e P 



v = vo 



M~ 2 



rd r ln(vo) a r 
( P . \ r 

IP 




d r P 



(15) 
(16) 



where we define the Mach number of the stationary solution M 2 o = ~ 10 3 (ijjj) ( " q+1 \ 



Thus, considering a polytropical gas, the only prescription of cr(r, 9) gives a quasi stationary solution. 
Our vortex model is given by: 



a — 1 + A exp 



UJ 2 



r 



2 

i) + 



i 

x 2 



(0 - Oof 



(17) 



where ro and #o define the vortex position, A the amplitude at the vortex center, u r the radial width 
and x the aspect ratio of the vortex. 

3 Numerical simulations 

The temporal solving of the Euler equations is made with a numerical scheme using the finite volume 
method, an exact Riemman solver, and in which the calculation method of the radial flux was improved 
regarding the classical MUSCL scheme. The time integration is based on the 2d order Runge-Kutta 
method. The temporal stability is verified which allows to integrate the solution for more than 1000 
disk rotations without significant errors (10 -8 on the speeds and 10 -6 on the density). The boundary 
conditions are the continuity to the stationary solution ao, Po, Vo m the radial direction, and periodic 
conditions in the azimuthal direction. 

We use a polar grid, with n r — 500 cells in the radial directio and no — 1000 cells in the azimuthal 
direction. Moreover four ghost cells are used outside the domain to prevent any perturbation due to 
the boundaries. Finally, the disk extends from 5 AU to 10 AU. 

3.1 The Rossby Waves Instability 

The first set of simulations is the study of the so-called Rossby Waves Instability. We will not look 
closely into it so we refer to the articles [8] and [9j. Those studies emphasize the fact that an axisy- 
metrical perturbation with a Gaussian radial profile can trigger it if the height and the width of this 
perturbation are sufficient. Hence we introduce the following initial conditions in our simulations: 

(18) 
(19) 
(20) 

(21) 
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Pini 


= (i + /)' 


^ini 
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= M- 2 
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where / is the gaussian perturbation. We set f = fo exp 



(r-r ) 2 



performed 6 simulations whith the following parameters fo and ujf. 



i + / 

and ro = 7.5 AU. We have 
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(bl) (b2) (b3) 

Figure 1: Density and vorticity relatively to the stationary solution at the initial time (al et bl), 
after 16 rotations (a2 et b2) and after 180 rotations (a3 et b3) for the parameters fo = 20% and 
oj f = 0.3 AU 



3.2 Vortex models 

We have studied several vortex models with different parameters. The vortex position is ro = 7.5 AU 
et #o = 3 rad. Three types of simulations have been performed: 
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4 Discussion 

The evolution of vortices is similar either they are issued from the Rossby Waves Instability or imposed 
in the disk. First, the several simulations confirm that vortices are robust solutions of the Euler 
equations. The Figures [T] and [2j which show the time evolution of different simulations, highlight that 
these structures conserve their mass and their vorticity for more than a hundreds rotations. They 
continue to exist during more than 10 4 years, which represents an important fraction of the disk life. 
As a consequence, it can not be considered as transitory solutions. The range of solutions is vaste in 
terms of radial width and of aspect ratios. The Figure [3] shows these different parameters concerning 
different simulations. Only the radial width is more constrained. Indeed the shear of the surrounding 
gas restricts the radial width to the disk scale height: 

TC 

H = — — ~ 0.45 AU a 7.5 AU (22) 

Nevertheless, most of the vortices have a radial width uj r wider than 2H and are stable ([7]). It seems 
that the flow is modified when the vortex is present, thus the structure undergoes a weak shearing 
between its inner and outer boundaries. The most noticeable phenomenon is the radial migration 
of vortices. They loose their kinetic momentum with a constant rate as the Figure [4] illustrates. 
The emission of density waves whose spiral pattern is due to the Keplerian shear, is systematic and 
responsible for the kinetic momentum transport. The article [10] proposes a linear study of the waves 
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Figure 2: Density and vorticity relatively to the stationary solution at the initial time (al et bl) and 
after 200 rotations (a2 et b2) for a model with the parameters A — 60%, uj r — 0.6 AU and x — 10 





(a) 



(b) 



Figure 3: Time evolution of u r (a) and \ (b) for different simulations 
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(a) (b) 
Figure 4: Time evolution of the vortex position (a) and the vortex mass (b) for different simulations 




Radius (AU) Wave action -s 



(a) (b) 

Figure 5: Radial transport of kinnetic momentum (a) and its variation with the migration velocity 
(b) for different simulations 

emitted by a perturbation and introduces the term of radial flux of the kinetic momentum which is 
written like the following equation: 

F r = 27TT < ravu > (23) 

where < > stands for the azimuthal average. Such a flux is not linked with linear waves. Nevertheless, 
we emphasize that the waves of our simulated vortices have radial transport which is constant in time. 
The wave of the inner part has an inferior rate than the one of the outer wave (Figure [5}a). The 
vortex situated at the center of the gap loses kinetic momentum. This radial flux is directed to the 
outside, as a consequence the disk and the vortex migrate onto the star. We illustrate Figure [5}b that 
this transport inside the waves is proportional to the vortex migration speed. 

The vortices migration speed is situated between 0.2 — 2 — . This migration appears thus to be 
fundamental in planet formation models and disk evolution scenarios. This migration can justify hot 
jupiters formation or a model in which the accretion is kept going by a huge number of vortices split 
in the disk. The vortex model that we have developed allows us to better understand this non-linear 
phenomenon and we plan to present the continuation of our work in a future publication. 
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